Van der Waals density functional: an appropriate exchange functional 
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In this paper, an exchange functional which is compatible with the non-local Rutgers-Chalmers 
correlation functional (vdW-DF) is presented. This functional, when employed with vdW-DF, 
demonstrates remarkable improvements on intermolecular separation distances while further im- 
proving the accuracy of vdW-DF interaction energies. The key to the success of this three param- 
eter functional is its reduction of short range exchange repulsion through matching to the gradient 
expansion approximation in the slowly varying/high density limit while recovering the large reduced 
gradient, s, limit set in the revised PBE exchange functional. This augmented exchange functional 
could be a solution to long-standing issues of vdW-DF lending to further applicability of density 
functional theory to the study of relatively large, dispersion bound (van der Waals) complexes. 
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Van der Waals, or London dispersion, interactions have 
profound importance in bio-organic systems as well as 
many novel materials being investigated for energy ap- 
plications. Despite the importance of these systems 
and applications, first principles simulations have been 
greatly lacking. The primary reason for this is the in- 
ability of traditional density functional theory exchange- 
correlation functionals to account for long-ranged, van 
der Waals interactions. This has limited first-principles 
investigations to quantum chemical methods which, due 
to their computational expense, are only capable of mod- 
eling fragments of the true material; thus often overlook- 
ing some of the more salient features of these systems. 

A promising solution to the problem of balancing speed 
and scalability with accuracy lies in the non-local cor- 
relation functional of the Rutgers-Chalmers collabora- 
tion, the aptly named van der Waals density functional 
(vdW-DF).0,"0| This method includes long-range disper- 
sion effects as a simple perturbation to the local den- 
sity approximation correlation term and has been ex- 
tremely successful in describing a diverse group of ma- 
terials properties - ranging from molecules to bulk poly- 
mers and the adsorption of molecules to surfaces and 
within bulk materials. Q Recent developments show that 
self-consistency gives no appreciable differences in com- 
puted interaction energies and such a non-local func- 
tional can in fact be incorporated in an extremely ef- 
ficient manner. [J] However, the overwhelming success of 
vdW-DF is marred by its consistent overestimation of in- 
termolecular distances. Q Analysis of various generalized 
gradient approximation exchange functionals (GGA X ) in- 
dicate that traditional functionals are either too repulsive 
at short distances or incorrectly exhibit some "correla- 
tion" binding at larger distances. The standard func- 
tional used within the vdW-DF, the revised Perdew- 
Burke-Erzenhoff functional (revPBE), [B| unfortunately 
gives too much repulsion at short distances. Replac- 
ing revPBE with Hartree-Fock (HF) exchange shows im- 
provements in the inter-species separation distances ob- 



tained with vdW-DF but at the cost of over-binding, i.e. 
considerably larger interaction energies than obtained via 
CCSD(T). [f| [7j Recent work suggests that for many dis- 
persion bound materials the PW86 functional [g[ most 
closely matches HF exchange, [il Hoj Similarly when ap- 
plied with vdW-DF it also strongly overbinds (see Fig. [2]). 

In this paper, a GGA X that may be more suitable for 
use with the vdW-DF correlation functional is proposed. 
This functional is derived through the introduction of an 
enhancement factor which obeys two specific constraints: 

(i) matching to the gradient expansion approximation 
(GEA) in the slowly varying/high density limit and 

(ii) a smooth asymptote to the upper bound empirically 
set in revPBE exchange. Initial results indicate dramatic 
improvements in vdW-DF separation distances while re- 
taining the accuracy of this method for a range of sys- 
tems. Most notable are improvements in the interaction 
energies and the intermolecular/interplanar se par ation 
distances obtained for S22 database structures [l2| and 
graphite. 

The general formula of a GGA X can be written as: 



E? GA = / d*rne x nll {n)F x (s), 



(1) 



where e x (n) is the exchange energy per particle in a 
uniform gas (e" nif (n) = -3ek F /iiT with k F = [i^n] 1 / 3 ) 
and F x (s) is the enhancement factor which is a function of 
s = \7n/(2kpn). This form of exchange ensures proper, 
uniform density scaling [l3| where F x (s) = 1 simply gives 
LDA exchange. In general, the enhancement factor is 
chosen such that F x (0) = 1. Here we design an F x (s) to 
fulfill two further criteria: 

(i) To reduce the short range exchange repulsion, in 
the limit of s — > 0, i.e. for slowly varying/high densities, 
the functional approaches the GEA: 
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FIG. 1: (color online)Enhancement factor for various GGA^ 
functionals. The blue, red and green dashed lines represent 
the enhancement factor of PBE, revPBE Q and PW86,[|| 
respectively. The black, dashed-dot line is the GEA and the 
orange, open circles are points from the proposed C09 x . The 
solid, black line indicates the revPBE upper bound. 
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FIG. 2: (color online) Interaction energy as a function of sepa- 
ration distance for the benzene dimer in the stacked sandwich 
configuration. The vdW-DF results with revPBE, C09 and 
PW86 for exchange are represented by a black dashed line 
with open circles, a red solid line with closed circles and a 
solid blue line with closed triangles, respectively. Benchmark 
CCSD(T) data [TJ] (blue cross) are plotted for reference. 



where /i = 0.0864. [llj Figure [T] depicts the enhancement 
factor of a number of GGA X . It can be clearly seen that 
many of these functionals deviate quite quickly from the 
GEA. Note that this constraint is similar to that used 
in the recent PBEsol exchange functional which was de- 
signed to restore the gradient expansion in order to re- 
move artificial bias towards free atoms. [3] Here, decreas- 
ing the enhancement factor for small s (thus restoring the 
GEA) leads to a reduction in the short range repulsion 
in the GGA X . 

(ii) The second constraint used in the proposed GGA X 
is to asymptote the revPBE upper bound of 2.245 in the 
large s- limit. [H This simple bound is taken to be compat- 
ible with previous applications of the revPBE exchange 
functional with the non-local vdW-DF correlation term. 
Here, we find that an empirical F x (s) bound similar to 
revPBE gives the best interaction energies. 

Using these constraints a simple, smooth, function can 
be constructed in the form: 



F x (s) = 1 



/is 2 e 



k(1 
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(3) 



with ^=0.0617, k=1.245, a=0.0483. Fig. Q] displays the 
enhancement factor of eqn. [3] along with that for other 
GGAs. The parameters were determined by simultane- 
ously fitting Eqn. [3] to GEA for s < 1.5 and to revPBE 
for s > 8.0. This fitting domain was arbitrarily chosen 
to allow for a decrease in F x (s) for small s and a smooth 
recovery of revPBE for large values of s. (In accordance 
with previous naming conventions this functional shall 
be referred to as 009^). The complimentary exchange 
potential can be easily constructed from the functional 
derivative of Eqn. [3] as shown in Eqn. 24 of Ref. d. 

To test the compatibility of the proposed C09 x with 
vdW-DF, self-consistent calculations within a modified 



version of the Abinit plane wave code 16| were per- 
formed. All calculations were carried out with a 30 Ha 
planewave cutoff and a single k-point at T. To reduce 
the effects of periodic images, all simulation cells were 
padded with at least 10 A of vacuum in all directions. 

The interaction energy, AE mt , as a function of sepa- 
ration distance, d scp , for the benzene dimer stacked in 
the sandwich configuration is plotted in Fig. [2j A com- 
parison of vdW-DF with the standard revPBE exchange 
functional (vdW-DF rcvPBE ) and the exchange functional 
(vdW-DF C09x ) using the enhancement factor of eqn. [3] 
shows a significant shortening of the separation distance 
from 4.07 A to 3.87 A. The vdW-DF C09 - is now in much 
better agreement, with regards to both interaction en- 
ergy and separation distance, with both the benchmark 
CCSD(T) (AE int = 1.81 kcal/mol and d scp = 3.70 A) 
as well as SAPT(DFT) (AE int = 1.67 kcal/mol and d scp 
= 3.80 A) [24| calculations. As previously mentioned, 
the PW86 exchange functional, which was recently re- 
ported to mimic Hartree-Fock exchange for dispersion 
bound complexes, gives excellent separation distances, 
but significantly over estimates the interaction energy. 

To further illustrate the value of the C09 x , the inter- 
action energy as a function of c lattice parameter for 
graphite is plotted in Fig. [3l Here, C09 x offers signifi- 
cant improvements in both the value of graphite lattice 
constants as well the interplanar interaction energy. Our 
computed lattice constant of 6.56 A is within 2% of ex- 
periment (6.70 A).[2oT| The interplanar binding energy 
(-59 kcal/mo)is also greatly improved (-52 ± 5 kcal/mol 



for experiment 25]). Note, vdW-DF 



rcvPBE x 



gives inter- 



action energies and lattice constants of -39 kcal/mol and 
7.35 A, respectively. 

A more stringent evaluation of the accuracy of the 
functional can be gained through comparison with the 
benchmark S22 database of Jurecka et al.[12] This 
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TABLE I: Computed interaction energies for the S22 model data set.pjj] C09 values in parentheses are for full geometry 
optimizations. For comparison, vdW-DF interactions energies using the revPBE and PBE exchange functionals taken from 
Ref. Il3 (unless otherwise noted) are listed. Deformation energies are not included. Energies are in kcal/mol. 



No. 



Complex 



vdW-DF 
revPBE a PBE a C09 



Benchmark 
CCSD(T)/CBS 



Hydrogen bonded complexes 



1 (NH 3 ) 2 (C* 2h ) 


2.44 


3.71 


2.88 


(2.92) 


3.17 


2 (H 2 0) 2 (C s ) 


4.08 


5.58 


4.99 


(5.15) 


5.02 


3 Formic acid dimer (C 2 h) 


14.07 


18.28 


20.12 


(24.59) 


18.61 


4 Formamide dimer (C 2 h) 


12.50 


16.10 


16.15 


(17.15) 


15.96 


5 Uracil dimer (C 2 h) 


16.17 


20.59 


20.93 


(22.36) 


20.65 


6 2-pyridoxine ■ 2-aminopyridine (Ci) 


14.02 


17.94 


17.67 


(19.66) 


16.71 


7 Adenine • thymine WC (Ci ) 


15. 19 b 


17.57 


17.25 


(19.16) 


16.37 


Complexes with predominant dispersion contribution 












8 (CH 4 ) 2 (£>3d) 


0.88 


1.55 


0.51 


(0.51) 


0.53 


9 (C 2 H 4 ) 2 (D 3d ) 


1.41 


2.68 


1.16 


(1.16) 


1.51 


10 Benzene • CH 4 (C 3 ) 


1.57 c 


2.51 


1.50 


(1.71) 


1.50 


11 Benzene dimer (C 2 h) 


2.74 d 


4.96 


3.32 


(3.40) 


2.73 


12 Pyrazine dimer (C s ) 


3.87 


6.25 


4.74 


(4.75) 


4.42 


13 Uracil dimer (C 2 ) 


9.4U 


12.91 


10.31 


(10.41) 


10.12 


14 Indole • benzene (Ci) 


A 'iA 
4.o4 


D.ZO 


Z.AA 
0.44 


( F. AQ\ 


K OO 


15 Adenine • thymine stack (Ci) 


10.60 e 


14.74 


12.79 


(12.73) 


12.23 


Mixed complexes 












16 Ethene • ethine (C 2v )) 


1.55 


2.38 


1.61 


(1.60) 


1.53 


17 Benzene • H 2 (C s ) 


2.72 f 


4.15 


3.25 


(3.24) 


3.28 


18 Benzene • NH 3 (C s ) 


1.87 


3.18 


2.28 


(2.27) 


2.35 


19 Benzene • HCN (C s ) 


3.87 


5.488 


4.51 


(4.48) 


4.46 


20 Benzene dimer (C 2v )) 


2.05 d 


3.98 


2.85 


(2.84) 


2.74 


21 Indole • benzene T-shape (Ci) 


4.72 


6.90 


5.75 


(5.71) 


5.73 


22 Phenol dimer (Ci) 


5.81 


8.51 


7.00 


(7.20) 


7.05 


Avg. % deviation 


18 


36 


5 


(9) 





structures optimized for separation distance 



*> Ref. 
c Ref. 
d Ref.ra 
e Ref. 13 
•fRef.[H 



database contains the interaction energies and structures 
of 22 structures with varying degrees of hydrogen bond- 
ing and vdW interactions computed with CCSD(T) ex- 
trapolated to the complete basis set limit and is currently 
accepted as the gold standard for theoretical methods 
used to study systems with significant dispersion inter- 
actions. Recently, Gulans and coworkers examined the 
S22 database using vdW-DF with both the revPBE and 
PBE exchange functionals. [li| In general, they found rea- 
sonable agreement with the interaction energies of the 
S22 database; with PBE producing better results for hy- 
drogen bonded complexes and revPBE showing less de- 
viation for dispersion dominated interactions. However, 
there results were all for vdW-DF optimized structures; 
which always give too large separation distances. 



Table U lists the computed interaction energies for 
the S22 database. Here, it is evident that across the 
board vdW-DF C09x is in much better agreement with 
the CCSD(T) benchmark values than when vdW-DF is 
used with either the PBE or revPBE functionals. In 
fact, vdW-DF C09x has an average percent deviation of 5% 
(9% if full geometry optimizations were performed); far 
less than revPBE and PBE (17% and 36%, respectively). 
Even more important is the fact that these interaction 
energies were obtained using the published S22 geome- 
tries - without any adjustment of dimer separation dis- 
tances, demonstrating once again the improvement that 
this functional offers with regards to both interaction en- 
ergies and determining optimum separation distances. 

All 22 structures were subsequently relaxed such that 
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FIG. 3: (color online) Interaction energy as a function of inter- 
planar separation distance for graphite. The vdW-DF results 
with revPBE and C09 for exchange are represented by a black 
dashed line with open circles and a solid black line with closed 
triangles, respectively. A recent experimental value [25[ (blue 
open circles with error bars) is plotted for reference. 



the forces on all the atoms were less than 0.02 eV/A. 
These values are listed in parentheses in the Table |U 
Analysis of the relaxed geometries indicate that the ma- 
jority of the deviations are related to changes in internal 
bond lengths. This is evident in the larger changes in the 
vdW-DF C09x interaction energies for hydrogen-bonded 
structures. 

In summary, an exchange functional that is compat- 
ible with the Rutgers-Chalmers van der Waals correla- 
tion functional is proposed. This functional was derived 
to closely match the enhancement factor F(s) of the 
gradient exchange approximation for values of < s < 
1.5, while having an asymptote to the revPBE bound of 
2.245. In general, this functional shows significant im- 
provements over the previous revPBE exchange. In par- 
ticular, vdW-DF C09x offers better agreement with the 
benchmark S22 database with an average deviation of 
only 5% at the intermolecular separation distances of 
the published geometries. This is a feat which far sur- 
passes that of previous vdW-DF calculations which re- 
quired larger separation distances. It should be noted 
that as previously pointed out that while the restoration 
of the gradient expansion gives improved dispersion inter- 
actions it is expected to worsen atomization energies. [l5| 
Nevertheless, these results highlight the promise of this 
functional for use with the vdW-DF method and may of- 
fer a pathway to even more accurate first-principles cal- 
culations of dispersion bound systems. 
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